$TITLE  GTAPinGAMS -- Static Multiregional Core Model in Algebraic MCP Syntax

*	Note:

*		This is the model implemented as a nonlinear program.

*		This implementation accomodates both constant-elasticity of
*		transformation between production for domestic and export
*		markets (eta < +INF), and perfect substitution between 
*		those markets (eta=+INF).

*		Variables, equations and GAMS keywords are in UPPER case.
*		Sets and parameters are in lower case.

*	Read the dataset:

$if not set dataset $set dataset gtap5_small
$include mrtdata

SCALAR
	eta	Elasticity of transformation - domestic vs. exports	/ +inf /,
	esubdm	Elasticity of substitution - domestic vs. imports	/ 4 /,
	esubmm  Elasticity of substitution - imports			/ 8 /;

ALIAS (f,ff);

PARAMETER
	vad(i,r)	Sectoral value-added,
	tau(i,r,s)	Unit transport cost coefficient

	thetaf(f,i,r)	Value added,
	thetad(i,r)	Domestic output,
	thetag(i,r)	Government demand, 
	thetap(i,r)	Private demand, 
	thetat(i,r)	Transport, 
	thetam(d,i,r)	Import value share,
	beta(i,s,r)	Value share of bilateral imports,
	gamma(i,s,r)	Goods share of unit import cost;

vad(i,r) = sum(f, vfm(f,i,r)*pf0(f,i,r));
tau(i,r,s)$vxmd(i,r,s) = vtwr(i,r,s) / vxmd(i,r,s);

thetam(d,i,r)$va(d,i,r) = vm(d,i,r) / va(d,i,r);
thetaf(f,i,r)$vad(i,r) = pf0(f,i,r) * vfm(f,i,r) / vad(i,r);
thetad(i,r)$vom(i,r) = vdm(i,r) / vom(i,r);
thetag(i,r) = pg0(i,r) * vgm(i,r) / vg(r);
thetap(i,r) = pc0(i,r) * vpm(i,r) / vp(r);
thetat(i,r)$sum((j,s), vst(j,s)) = vst(i,r) / vt;
beta(i,s,r)$vxmd(i,s,r) = (vxmd(i,s,r)*pmx0(i,s,r)+vtwr(i,s,r)*pmt0(i,s,r) ) / vim(i,r);
gamma(i,s,r)$vxmd(i,s,r) = vxmd(i,s,r)*pmx0(i,s,r) /
			  (vxmd(i,s,r)*pmx0(i,s,r) + vtwr(i,s,r)*pmt0(i,s,r) );

*	Declare subsets of d in order to make it possible to write
*	equation MKT_PA:

SETS	c_d(d)/c/,  g_d(d)/g/, i_d(d) /i/;

VARIABLES
	OBJ		Dummy objective function

	A_G(i,r)	Public sector unit demand
	A_C(i,r)	Private unit demand
	A_F(f,i,r)	Factor unit demand
	A_X(i,r)	Export unit supply
	A_D(i,r)	Domestic unit supply
	A_M(i,r,s)	Import unit demand

	C(r)		Private consumption
	G(r)		Public provision
	Y(i,r)		Aggregate Output
	M(i,r)		Import aggregation
	A(D,i,r)	Armington aggregation of domestic and imports
	YT		Transport 

	PC(r)		Private demand
	PG(r)		Public provision
	PD(i,r)		Domestic price
	PX(i,r)		Export price
	PY(i,r)		Output price (for ETA=INF)
	PM(i,r)		Import price
	PA(D,i,r)	Armington composite price
	PF(f,r)		Factor price
	PT		Transport services

	RA(r)		Representative agent income;

EQUATIONS

	OBJDEF		Defines OBJ,

	DEF_G(i,r)	Public sector demand
	DEF_C(i,r)	Private demand
	DEF_F(f,i,r)	Factor demand
	DEF_X(i,r)	Export supply
	DEF_D(i,r)	Domestic supply
	DEF_M(i,r,s)	Import demand

	PRF_C(r)	Private consumption
	PRF_G(r)	Public provision
	PRF_Y(i,r)	Aggregate Output
	PRF_M(i,r)	Import aggregation
	PRF_A(D,i,r)	Armington aggregation of domestic and imports
	PRF_YT		Transport 

	MKT_PC(r)	Private demand
	MKT_PG(r)	Public provision
	MKT_PY(i,r)	Output price
	MKT_PD(i,r)	Domestic price
	MKT_PX(i,r)	Export price
	MKT_PM(i,r)	Import price
	MKT_PA(D,i,r)	Armington composite price
	MKT_PF(f,r)	Factor price
	MKT_PT		Transport services

	INC_RA(r)	Representative agent;


*	Compensated unit demand and supply functions:
*	============================================

DEF_G(i,r)$vgm(i,r)..

	A_G(i,r) =E= vgm(i,r) * 
	  PROD(j, (PA("g",j,r)*(1+tg(j,r))/pg0(j,r))**thetag(j,r)) /
	  ( PA("g",i,r)*(1+tg(i,r))/pg0(i,r) );

DEF_C(i,r)$vpm(i,r)..

	A_C(i,r) =E= vpm(i,r) * 
	  PROD(j, (PA("c",j,r)*(1+tp(j,r))/pc0(j,r))**thetap(j,r)) /
	  ( PA("c",i,r)*(1+tp(i,r))/pc0(i,r) );

DEF_F(f,i,r)$vfm(f,i,r)..

	A_F(f,i,r) =E= vfm(f,i,r) *
	  PROD(ff, (PF(ff,r)*(1+tf(ff,i,r))/pf0(ff,i,r) )**thetaf(ff,i,r)) /
	  ( PF(f,r)*(1+tf(f,i,r)) / pf0(f,i,r) );

DEF_X(i,r)$(vxm(i,r) and (1/eta gt 0)) ..

	A_X(i,r) =E=  vxm(i,r) * (PX(i,r) / 
	    (thetad(i,r)  * PD(i,r)**(1+eta) + 
	  (1-thetad(i,r)) * PX(i,r)**(1+eta))**(1/(1+eta)) )**eta;

DEF_D(i,r)$(vdm(i,r) and (1/eta gt 0))..

	A_D(i,r) =E=  vdm(i,r) * (PD(i,r) / 
	    (thetad(i,r)  * PD(i,r)**(1+eta) + 
	  (1-thetad(i,r)) * PX(i,r)**(1+eta))**(1/(1+eta)) )**eta;

DEF_M(i,r,s)$vxmd(i,r,s)..

	A_M(i,r,s) =E= vxmd(i,r,s) * ( PM(i,s) / 
	     ( gamma(i,r,s)  * (PX(i,r)$(1/eta gt 0) + PY(i,r)$(1/eta eq 0))
		*(1+tx(i,r,s))*(1+tm(i,r,s))/pmx0(i,r,s)
	   + (1-gamma(i,r,s)) * PT*(1+tm(i,r,s)) / pmt0(i,r,s) ) )**esubmm;


*	Zero profit per unit activity are defined here:
*	==============================================

*	Production:

PRF_Y(i,r)$vom(i,r)..

	SUM(j, vafm(j,i,r) * PA("i",j,r) * (1+ti(j,i,r)) ) +
	SUM(f, A_F(f,i,r) * PF(f,r) * (1 + tf(f,i,r)) )
        =E= (1 - ty(i,r)) * ((PY(i,r) * vom(i,r))$(1/eta eq 0)
		   + (PD(i,r) * A_D(i,r) + PX(i,r) * A_X(i,r))$(1/eta gt 0));

*	Import-domestic aggregation by submarket:

PRF_A(d,i,r)$va(d,i,r)..

	  ((1-thetam(d,i,r)) * (PD(i,r)$(1/eta) + PY(i,r)$(1/eta eq 0))**(1-esubdm) +
	      thetam(d,i,r)  * PM(i,r)**(1-esubdm) )**(1/(1-esubdm)) =E= PA(d,i,r);

*	Armington aggregation across imports from different countries:

PRF_M(i,r)$vim(i,r)..

	SUM(s, (1 + tm(i,s,r)) * A_M(i,s,r) * 
	((PX(i,s)$(1/eta gt 0)+PY(i,s)$(1/eta eq 0)) * (1 + tx(i,s,r)) 
		+ PT * tau(i,s,r))) =E= PM(i,r) * vim(i,r) ;

*	Public output:

PRF_G(r)..	SUM(i, PA("g",i,r) * (1+tg(i,r)) * A_G(i,r)) =E= PG(r) * vg(r);

*	Private consumption:

PRF_C(r)..	SUM(i, PA("c",i,r) * (1+tp(i,r)) * A_C(i,r)) =E= PC(r) * vp(r);


*	Inter-national transport services (Cobb-Douglas):

PRF_YT..	PROD((i,r), (PX(i,r)$(1/eta gt 0) + PY(i,r)$(1/eta eq 0))**thetat(i,r)) =E= PT;

*	Market clearing
*	===============

*	Aggregate output (when domestic and exports are perfect substitutes):

MKT_PY(i,r)$(vom(i,r) and (1/eta=0))..

	vom(i,r) * Y(i,r) =e= SUM(s, A_M(i,r,s) * M(i,s)) + VST(i,r) * YT * (PT/PY(i,r))
	  +  SUM(d, A(d,i,r) * vd(d,i,r) * ( PA(d,i,r)/PY(i,r) )**esubdm ) + vi(r)$cgd(i);


*	Exports:

MKT_PX(i,r)$(vxm(i,r) and (1/eta gt 0))..

	A_X(i,r) * Y(i,r) =E= SUM(s, A_M(i,r,s) * M(i,s)) + VST(i,r) * YT * (PT/PX(i,r));

*	Domestic supply:

MKT_PD(i,r)$(vdm(i,r) and (1/eta gt 0))..

	A_D(i,r) * Y(i,r) =E= SUM(d, A(d,i,r) * vd(d,i,r) * ( PA(d,i,r)/PD(i,r) )**esubdm )
				+ vi(r)$cgd(i);

*	Imports:
 
MKT_PM(i,r)$vim(i,r)..

	vim(i,r) * M(i,r) =E= 
		SUM(d, A(d,i,r) * vm(d,i,r) * ( PA(d,i,r)/PM(i,r) )**esubdm );		 
  
*	International transport:

MKT_PT..	YT * vt =E=  sum((i,r,s), A_M(i,r,s) * M(i,s) * tau(i,r,s));

*	Armington supply:

MKT_PA(d,i,r)$va(d,i,r)..

	va(d,i,r) * A(d,i,r) =E= sum(j, vafm(i,j,r) * Y(j,r))$i_d(d)

			+ (A_C(i,r) * C(r))$c_d(d) + ( A_G(i,r) *  G(r))$g_d(d);
		      
*	Government provision:

MKT_PG(r)..	G(r) =E= 1;

*	Factor market:

MKT_PF(f,r)..	evoa(f,r) =E= sum(i, A_F(f,i,r) * Y(i,r));

	
*	Private demand:

MKT_PC(r)..	C(r) * vp(r) =E= RA(r) / PC(r) ;

*	Income balance
*	==============

INC_RA(r)$(RA.LO(r) ne RA.UP(r))..	

	RA(r) =E= sum(f,   PF(f,r) * evoa(f,r)) 
		  + sum(num, PC(num) * vb(r))
		  - sum(cgd, (PD(cgd,r)$(1/eta) + PY(cgd,r)$(1/eta eq 0))  * vi(r))
		  - PG(r) * vg(r)

* Output tax:

	+    sum(i, ty(i,r) * 

	  ( (PX(i,r) * A_X(i,r) + PD(i,r) * A_D(i,r))$(1/eta gt 0)
			       + (PY(i,r) * vom(i,r))$(1/eta eq 0) ) * Y(i,r))

* Tax on intermediate demand:

	+    sum((i,j), ti(j,i,r) * PA("i",j,r) * vafm(j,i,r) * Y(i,r) )

* Taxes on factor use:

	+    sum( (i,f),  tf(f,i,r) * PF(f,r) * A_F(f,i,r) * Y(i,r))

* Export tax:

	+     sum((i,s), tx(i,r,s) * A_M(i,r,s) * M(i,s) 

		* (PX(i,r)$(1/eta gt 0) + PY(i,r)$(1/eta eq 0)) )

* Import tariff applies to mechandise gross of export tax
* and transport cost:

	+    sum((i,s), tm(i,s,r) * A_M(i,s,r) * M(i,r) * 
		(  (PX(i,s)$(1/eta gt 0) + PY(i,s)$(1/eta eq 0))
			* (1+tx(i,s,r)) + PT * tau(i,s,r) ) )
       
* Taxes on goverment consumption:

	+    sum(i, tg(i,r) * PA("g",i,r) * A_G(i,r) * G(r))
		
* Taxes on private consumption:		    

	+    sum(i, tp(i,r) * PA("c",i,r) * A_C(i,r) * C(r));

OBJDEF..	OBJ =E= 0;

*	Define the model:

MODEL GTAP /
	DEF_G, DEF_C, DEF_F, DEF_X, DEF_D, DEF_M,
	PRF_C, PRF_G, PRF_Y, PRF_M, PRF_A, PRF_YT, 
	MKT_PC, MKT_PG, MKT_PD, MKT_PX, MKT_PY,
	MKT_PM, MKT_PA, MKT_PF, MKT_PT, 
	INC_RA, OBJDEF /;	


* Install lower bounds to avoid bad function calls:

PC.LO(r)=1.e-5;
PG.LO(r)=1.e-5; 
PD.LO(i,r)=1.e-5;
PX.LO(i,r)=1.e-5;
PY.LO(i,r)=1.e-5;
PM.LO(i,r)=1.e-5;
PA.LO(d,i,r)=1.e-5;
PF.LO(f,r)=1.e-5; 
PT.LO=1.e-5;

* Install level values:

C.L(r)=1;
G.L(r)=1;
Y.L(i,r)=1;
M.L(i,r)=1;
A.L(d,i,r)=1;

YT.L=1;
PC.L(r)=1;
PG.L(r)=1;
PD.L(i,r)=1;
PX.L(i,r)=1;
PY.L(i,r)=1;

PM.L(i,r)=1;
PA.L(d,i,r)=1;
PF.L(f,r)=1;
PT.L=1;

RA.L(r)=vp(r);

*	Fix a numeraire to simplify comparison with MGE:

RA.FX(num) = RA.L(num);
	
A_G.FX(i,r)$(not vgm(i,r)) = 0;
A_C.FX(i,r)$(not vpm(i,r)) = 0;
A_F.FX(f,i,r)$(not vfm(f,i,r)) = 0;
A_X.FX(i,r)$(not vxm(i,r)) = 0;
A_D.FX(i,r)$(not vdm(i,r)) = 0;
A_M.FX(i,r,s)$(not vxmd(i,r,s)) = 0;

Y.FX(i,r)$(not vom(i,r)) = 0;
M.FX(i,r)$(not vim(i,r)) = 0;
A.FX(d,i,r)$(not va(d,i,r)) = 0;
PA.FX(d,i,r)$(not va(d,i,r)) = 0.0001;
PD.FX(i,r)$(not vdm(i,r)) = 0.0001;
PX.FX(i,r)$(not vxm(i,r)) = 0.0001;
PM.FX(i,r)$(not vim(i,r)) = 0.0001;
PY.FX(i,r)$(not vom(i,r)) = 0.0001;

A_G.L(i,r) = vgm(i,r);
A_C.L(i,r) = vpm(i,r);
A_F.L(f,i,r) = vfm(f,i,r);
A_X.L(i,r) =  vxm(i,r);
A_D.L(i,r) =  vdm(i,r);
A_M.L(i,r,s) = vxmd(i,r,s);

*	Treat fixed variables as constants:

GTAP.HOLDFIXED = 1;

*	Give enough memory:

GTAP.WORKSPACE = 20;

* Replicate the benchmark

GTAP.ITERLIM=8000;
SOLVE GTAP USING NLP MAXIMIZING OBJ;

